\section{Navier Stokes Equation in Orthogonal Curvilinear Coordinate System}

Problem comes from the idea that we want to solve invisid fluid governing equation on arbitary smooth manifold \textbf{M} which we assume it to be diffeomorphic with 3D Euclidean space. We can definitely map a small region of \textbf{M} onto a chart with parameterization $c_1, c_2, c_3$.  Suppose the Euclidean parameterization is $x, y, z$. Then locally we have a bijection between those two parameterization.

\begin{equation}\nonumber
\begin{split}
x = \hat{X}(c_1, c_2, c_3)\\
y = \hat{Y}(c_1, c_2, c_3)\\
z = \hat{Z}(c_1, c_2, c_3)\\
c_1 = \hat{C}_1(x, y, z)\\
c_2 = \hat{C}_2(x, y, z)\\
c_3 = \hat{C}_3(x, y, z)\\
\end{split}
\end{equation}

Notice $x, y, z$ is actually also the coordinate function in Euclidean space, so a point can be represented by a vector \textbf{r},

\begin{equation}\nonumber
\textbf{r} = x\textbf{i} + y\textbf{j} + z\textbf{k}
\end{equation}

As we know in differential geometry, for any point \textbf{r} in 3D Euclidean space $\textbf{R}^3$, there is a tangent space $\textbf{T}_{\boldsymbol{r}}\textbf{R}^{3}$ which is diffeomorphic with $\textbf{T}_{\boldsymbol{0}}\textbf{R}^{3}$ on the origin. Thus it is also a vector space and can be spanned by basis vectors. If we take derivative of Euclidean coordinate functions with respect to our new parameterization, we can get the new basis vectors of each new parameter respectively (represented by Euclidean basis).

\begin{equation}\nonumber
\begin{split}
\textbf{r}_1 = \frac{\partial \textbf{r}}{\partial c_1} = \frac{\partial x}{\partial c_1}\textbf{i} + \frac{\partial y}{\partial c_1}\textbf{j} + \frac{\partial z}{\partial c_1}\textbf{k}\\
\textbf{r}_2 = \frac{\partial \textbf{r}}{\partial c_2} = \frac{\partial x}{\partial c_2}\textbf{i} + \frac{\partial y}{\partial c_2}\textbf{j} + \frac{\partial z}{\partial c_2}\textbf{k}\\
\textbf{r}_3 = \frac{\partial \textbf{r}}{\partial c_3} = \frac{\partial x}{\partial c_3}\textbf{i} + \frac{\partial y}{\partial c_3}\textbf{j} + \frac{\partial z}{\partial c_3}\textbf{k}\\
\end{split}
\end{equation}

To make our life easier, we will assume the basis vectors are orthogonal and actually we can always stick to Frenet framework for arbitary surface. However, one can also show it for arbitary skew basis with more patience cooking things up. Notice $\textbf{r}_i$ is not normalized yet, so after normalization we get the real basis vectors.

\begin{equation}\nonumber
\begin{split}
\textbf{e}_i = \frac{\textbf{r}_i}{\|\textbf{r}_i\|} &= \frac{\textbf{r}_i}{h_i}, h_i = \|\textbf{r}_i\|\\
\textbf{e}_i \cdot \textbf{e}_j &= \delta_{ij}
\end{split}
\end{equation}

We also need second order derivative of $\textbf{r}$ later,

\begin{equation}\nonumber
\begin{split}
\frac{\partial \textbf{r}_i}{\partial c_j} &= \textbf{r}_{ij}\\
&= \frac{\partial (h_i\textbf{e}_i)}{\partial c_j}\\
&= \frac{\partial h_i}{\partial c_j} \textbf{e}_i + h_i\frac{\partial \textbf{e}_i}{\partial c_j}\\
&= h_{ij}\textbf{e}_i + h_i\textbf{e}_{ij}
\end{split}
\end{equation}

It's not hard to prove $\textbf{e}_i \perp \textbf{e}_{ij}$ like how we did in \textbf{BTN} Frenet framework. Moreover, one can think of $\textbf{r}_i$ as a velocity vector. Its derivative is acceleration which can be decomposed into tangential and normal components. Here $\textbf{e}_{ij}$ is the normal direction.

Then we can derive all the second order derivative of the basis vectors $\textbf{e}_{i}$. Here we just show two cases of them. Others can be derived following the same idea.

\begin{equation}\nonumber
\begin{split}
\textbf{r}_1\cdot\textbf{r}_1 &= h_1^2\\
\textbf{r}_1\cdot\textbf{r}_{11} &= h_1h_{11}\\
\textbf{r}_1\cdot\textbf{r}_{12} &= h_1h_{12}\\
\textbf{r}_1\cdot\textbf{r}_{13} &= h_1h_{13}\\
\textbf{r}_1\cdot\textbf{r}_{2} &= 0\\
\textbf{r}_1\cdot\textbf{r}_{3} &= 0
\end{split}
\end{equation}
\begin{equation}\nonumber
\begin{split}
&\Rightarrow \frac{\partial \textbf{r}_1\cdot\textbf{r}_2}{\partial c_1} = 0\\
&\Rightarrow \textbf{r}_{11}\cdot\textbf{r}_2 + \textbf{r}_1\cdot\textbf{r}_{21} = 0\\
&\Rightarrow \textbf{r}_{11}\cdot\textbf{r}_{2} = -\textbf{r}_1\cdot\textbf{r}_{12}\\
&\Rightarrow \textbf{r}_{11}\cdot\textbf{r}_2 = -h_1h_{12}\\
&\Rightarrow \textbf{r}_{11}\cdot\textbf{r}_3 = -h_1h_{13}\\
&\Rightarrow \textbf{r}_{11} = h_{11}\textbf{e}_1 - \frac{h_1h_{12}}{h_2}\textbf{e}_2 - \frac{h_1h_{13}}{h_3}\textbf{e}_3 = h_{11}\textbf{e}_1 + h_1\textbf{e}_{11}\\
&\Rightarrow \textbf{e}_{11} = - \frac{h_{12}}{h_2}\textbf{e}_2 - \frac{h_{13}}{h_3}\textbf{e}_3
\end{split}
\end{equation}

\begin{equation}\nonumber
\begin{split}
\frac{\partial (\textbf{r}_i\textbf{r}_j)}{\partial c_k} = 0\\
\textbf{r}_{ij} = \textbf{r}_{ji}
\end{split}
\end{equation}
\begin{equation}\nonumber
\begin{split}
&\Rightarrow \textbf{r}_{13}\cdot\textbf{r}_2 + \textbf{r}_1\textbf{r}_{23} = 0\\
&\Rightarrow \textbf{r}_{21}\cdot\textbf{r}_3 + \textbf{r}_2\textbf{r}_{31} = 0\\
&\Rightarrow \textbf{r}_{32}\cdot\textbf{r}_1 + \textbf{r}_3\textbf{r}_{12} = 0\\
&\Rightarrow \textbf{r}_{32}\cdot\textbf{r}_1 = \textbf{r}_2\cdot\textbf{r}_{31}\\
&\Rightarrow \textbf{r}_{32}\cdot\textbf{r}_1 = - \textbf{r}_1\cdot\textbf{r}_{23}\\
&\Rightarrow \textbf{r}_{32}\cdot\textbf{r}_1 =0\\
&\Rightarrow \textbf{r}_{32} = h_{32}\textbf{e}_3 + h_3\textbf{e}_{32} = \textbf{r}_{23} = h_{23}\textbf{e}_2 + h_2\textbf{e}_{23}
\end{split}
\end{equation}

\begin{equation}\nonumber
\begin{split}
&\Rightarrow h_3\textbf{e}_{32} = h_{23}\textbf{e}_2\\
&\Rightarrow h_2\textbf{e}_{23} = h_{32}\textbf{e}_3\\
&\Rightarrow \textbf{e}_{32} = \frac{h_{23}}{h_3}\textbf{e}_2\\
&\Rightarrow \textbf{e}_{23} = \frac{h_{32}}{h_2}\textbf{e}_3
\end{split}
\end{equation}

Another problem to state before we get to our goal is the reprepsentation of gradient operator in the new frame, which is

\begin{equation}\nonumber
\nabla = \frac{1}{h_1}\frac{\partial}{\partial c_1}\textbf{e}_1 + \frac{1}{h_2}\frac{\partial}{\partial c_2}\textbf{e}_2 + \frac{1}{h_3}\frac{\partial}{\partial c_3}\textbf{e}_3
\end{equation}

One can certainly derive it from the invariant intrinsic property of the manifold. Another way to look at it is that gradient is covector like the basis vector. One can simply use the same O(3) group tensor which is applied to basis vectors and apply it to gradient vector.

\subsubsection*{Momentum Equation}

Now I start deriving the invisid Navier Stokes equation in the new manifold. The velocity field in the new manifold can be represented as 

\begin{equation}\nonumber
\textbf{u} = u_1\textbf{e}_1 + u_2\textbf{e}_2 + u_3\textbf{e}_3
\end{equation}

And the invisid N-S equation in Euclidean space is

\begin{equation}\nonumber
\frac{\partial \textbf{u}}{\partial t} + \textbf{u}\cdot\nabla\textbf{u} = -\frac{1}{\rho}\nabla p + \textbf{f}
\end{equation}

I will just expand the equation term by term.

$\bullet$ \textit{Time Derivative Term}

\begin{equation}\nonumber
\frac{\partial \textbf{u}}{\partial t} = \frac{\partial u_1}{\partial t}\textbf{e}_1 + \frac{\partial u_2}{\partial t}\textbf{e}_2 + \frac{\partial u_3}{\partial t}\textbf{e}_3
\end{equation}

$\bullet$ \textit{Advection Term}

Having the gradient operator on the new manifold, we can get the directional derivative operator directly

\begin{equation}\nonumber
\textbf{u}\cdot\nabla = \frac{u_1}{h_1}\frac{\partial}{\partial c_1} + \frac{u_2}{h_2}\frac{\partial}{\partial c_2} + \frac{u_3}{h_3}\frac{\partial}{\partial c_3}
\end{equation}

Then the advection becomes to

\begin{equation}\nonumber
\begin{split}
\textbf{u}\cdot\nabla\textbf{u} &= \textbf{u}\cdot\nabla(u_1\textbf{e}_1) + \textbf{u}\cdot\nabla(u_2\textbf{e}_2) + \textbf{u}\cdot\nabla(u_3\textbf{e}_3) \\
\textbf{u}\cdot\nabla(u_1\textbf{e}_1) &= \frac{u_1}{h_1}\frac{\partial (u_1\textbf{e}_1)}{\partial c_1} + \frac{u_2}{h_2}\frac{\partial (u_1\textbf{e}_1)}{\partial c_2} + \frac{u_3}{h_3}\frac{\partial (u_1\textbf{e}_1)}{\partial c_3} \\
&= \frac{u_1}{h_1}\frac{\partial u_1}{\partial c_1}\textbf{e}_1 + \frac{u_1^2}{h_1}\frac{\partial \textbf{e}_1}{\partial c_1} +       \frac{u_2}{h_2}\frac{\partial u_1}{\partial c_2}\textbf{e}_1 + \frac{u_1u_2}{h_2}\frac{\partial \textbf{e}_1}{\partial c_2} +      \frac{u_3}{h_3}\frac{\partial u_1}{\partial c_3}\textbf{e}_1 + \frac{u_1u_3}{h_3}\frac{\partial \textbf{e}_1}{\partial c_3} \\
&= (\frac{u_1}{h_1}\frac{\partial u_1}{\partial c_1} + \frac{u_2}{h_2}\frac{\partial u_1}{\partial c_2} + \frac{u_3}{h_3}\frac{\partial u_1}{\partial c_3})\textbf{e}_1 + (\frac{u_1^2}{h_1}\textbf{e}_{11} +  \frac{u_1u_2}{h_2}\textbf{e}_{12} + \frac{u_1u_3}{h_3}\textbf{e}_{13})\\
\textbf{u}\cdot\nabla(u_2\textbf{e}_2) &=  (\frac{u_1}{h_1}\frac{\partial u_2}{\partial c_1} + \frac{u_2}{h_2}\frac{\partial u_2}{\partial c_2} + \frac{u_3}{h_3}\frac{\partial u_2}{\partial c_3})\textbf{e}_2 + (\frac{u_1u_2}{h_1}\textbf{e}_{21} +  \frac{u_2^2}{h_2}\textbf{e}_{22} + \frac{u_2u_3}{h_3}\textbf{e}_{23})\\
\textbf{u}\cdot\nabla(u_3\textbf{e}_3) &=  (\frac{u_1}{h_1}\frac{\partial u_3}{\partial c_1} + \frac{u_2}{h_2}\frac{\partial u_3}{\partial c_2} + \frac{u_3}{h_3}\frac{\partial u_3}{\partial c_3})\textbf{e}_3 + (\frac{u_1u_3}{h_1}\textbf{e}_{31} +  \frac{u_2u_3}{h_2}\textbf{e}_{32} + \frac{u_3^2}{h_3}\textbf{e}_{33})
\end{split}
\end{equation}

By plugging the formula of $\textbf{e}_{ij}$, we can get to the advection term on the new manifold.

$\bullet$ \textit{Pressure Gradient Term}

\begin{equation}\nonumber
\nabla p = \frac{1}{h_1}\frac{\partial p}{\partial c_1}\textbf{e}_1 + \frac{1}{h_2}\frac{\partial p}{\partial c_2}\textbf{e}_2 + \frac{1}{h_3}\frac{\partial p}{\partial c_3}\textbf{e}_3
\end{equation}

$\bullet$ \textit{External Force Term}

\begin{equation}\nonumber
\textbf{f} = f_1\textbf{e}_1 + f_2\textbf{e}_2 + f_3\textbf{e}_3
\end{equation}